library(ggplot2)
library(tidyr)
library(dplyr)
library(patchwork)
library(ggridges)
library(devtools)
library(tidyverse)
library(wesanderson)
library(modelsummary)
library(ggpubr)

# leadership generation
opamed$leader <- NA
opamed$leader[opamed$leadership_transition == 0] <- 0
opamed$leader[opamed$leadership_transition > 0 & opamed$leadership_transition < 3] <- 1
opamed$leader[opamed$leadership_transition >= 3 & opamed$leadership_transition < 6] <- 2
opamed$leader[opamed$leadership_transition >= 6 ] <- 3
# decades
opamed$decade <- NA
opamed$decade[opamed$year > 1910 & opamed$year < 1950] <- 0
opamed$decade[opamed$year >= 1950 & opamed$year < 1960] <- 1
opamed$decade[opamed$year >= 1960 & opamed$year < 1970] <- 2
opamed$decade[opamed$year >= 1970 & opamed$year < 1980] <- 3
opamed$decade[opamed$year >= 1980 & opamed$year < 1990] <- 4
opamed$decade[opamed$year >= 1990 & opamed$year < 2000] <- 5
opamed$decade[opamed$year >= 2000 & opamed$year < 2010] <- 6
opamed$decade[opamed$year >= 2010] <- 7
# party history
opamed$party.history <- NA
opamed$party.history[opamed$foundation_miller < 26] <- 0
opamed$party.history[opamed$foundation_miller >= 26 & opamed$foundation_miller < 50] <- 1
opamed$party.history[opamed$foundation_miller >= 50] <- 2
# other computations
opamed$postwar <- ifelse(opamed$decade == 0, NA, opamed$decade)
opamed$dom.sing <- ifelse(opamed$regime < 2, opamed$regime, NA)
opamed$leaderDum <- ifelse(opamed$leader == 0, 0, 1)
opamed$gini <- (opamed$gini_swiid/100)
opamed$lnpop <- log(opamed$population_mln)
opamed$strong.party <- ifelse(opamed$opamed_strategy == 0, 1, 0)
opamed$lngini <- log(opamed$gini_swiid)
opamed$regime.dummy <- ifelse(opamed$regime == 0, 0, 1)
opamed$lnpop <- log(opamed$population_mln)
opamed$lngdp <- log(opamed$gdp_pc_mpd)
opamed$lnoil <- log(opamed$oil_pc)
opamed$lnoil <- ifelse(opamed$lnoil < -7.5, NA, opamed$lnoil)

## uploading MAPP dataset 
MAPP<- read_excel("Downloads/MAPP_dataset_-_Version_2.0.xlsx")
MAPP$membershipK <- (MAPP$`Number of party members (MAPP)`/1000)

# OPAMED summary table ####
data.summary <- opamed %>%
  select('Ruling Party' = party,
         'State' = state,
         'Party Membership Volatility' = memb_growth,
         'Ruling Party Size' = mempop_100,
         "Absolute Membership" = membership_1000,
         'Years Since Foundation' = foundation_miller,
         'Population (in millions)' = population_mln,
         "Log Population" = lnpop,
         'Survival Length' = Survival,
         "Party Survived (3-level)" = survived,
         "Elite Size WhoGov" = whogov_core,
         "Leadership Succession" = leadership_transition,
          "Election" = election,
         "Regime Type" = regime,
         "Autocratic Level" = p5_polity2,
         "Marxist" = marxist,
         'GDP per capita' = gdp_pc_mpd,
         "Log GDP per capita" = lngdp,
         'Oil per capita' = oil_pc,
         'Inequality Index SWIID' = gini_swiid,
         "GDP growth" = GDP_growth,
         "Year" = year)

datasummary_skim(data.summary, histogram=T)
# Figure 1. UMNO case study ####
# Create the data frame from OPAMED & electoral results
umno <- data.frame(
  election = c(2004,2008,2013, 2018),
  vote_share = c(2489 , 2382, 3252, 2526),
  membership = c(2400,2700, 3413, 3500), 
  stringsAsFactors = FALSE)

df <- umno %>%
  select(election, vote_share, membership) %>%
  gather(key = "variable", value = "value", -election)
head(df)

# color 
umno.plot <- ggplot(df, aes(x = election, y = value)) + 
  geom_line(aes(color = variable)) + 
  scale_color_manual(values = wes_palette("BottleRocket2", n = 2),
                     labels = c("Absolute membership", "Vote-share"))+
  scale_x_continuous(breaks=seq(2004, 2018, by = 1)) +
  xlab("Years") + ylab("Thousands of citizens") +
  labs(color = "UMNO") +
  theme_light()+
  theme(legend.position = "bottom")

tiff("Figure 1.tiff", units = "in", width = 8, height = 5, res = 350)
umno.plot
dev.off()

cairo_ps("Figure 1.eps")
umno.plot
dev.off()

# grey
umno.plot.grey <- ggplot(df, aes(x = election, y = value, linetype = factor(variable))) + 
  geom_line(aes(group = factor(variable))) +
  scale_x_continuous(breaks=seq(2004, 2018, by = 1)) +
  xlab("Years") + ylab("Thousands of citizens") +
  labs(linetype = "UMNO") +
  theme_light()+
  theme(legend.position = "bottom")

tiff("Figure 1 b&w.tiff", units = "in", width = 8, height = 5, res = 350)
  umno.plot.grey
dev.off()

cairo_ps("Figure 1 bw.eps")
umno.plot.grey
dev.off()

# Figure 2. MAPP & OPAMED Comparison ####
# A general plot of party memberhsip growth ###

opamed.plot <- ggplot(opamed, aes(x = year, y = mempop_100)) + 
  geom_point() +
  geom_smooth(method = "loess", se=T) +
  xlab("Years") + ylab("OPAMED ruling party penetration (% of population)") +
  theme_light() +
  theme(legend.position="bottom") 

## MAPP plot ##
  
mapp <- ggplot(MAPP, aes(x = `Year membership`, y = membershipK)) + 
  geom_point() +
  geom_smooth(method = "loess", se=T) +
  xlab("Years") + ylab("MAPP political party membership (in thousands)") +
  theme_light() +
  theme(legend.position="bottom") 

comparison <- mapp + opamed.plot
tiff("Figure 2.tiff", units = "in", width = 10, height = 10, res = 350)
comparison
dev.off()

cairo_ps("Figure 2.eps")
comparison
dev.off()

# Fig. 2 black and white
opamed.plotbw <- ggplot(opamed, aes(x = year, y = mempop_100)) + 
  geom_point(alpha = .3) +
  geom_smooth(method = "loess", se=T, color = "black", linetype = "dashed") +
  xlab("Years") + ylab("OPAMED ruling party penetration (% of population)") +
  theme_light() +
  theme(legend.position="bottom") 

## MAPP plot ##

mappbw <- ggplot(MAPP, aes(x = `Year membership`, y = membershipK)) + 
  geom_point(alpha = .3) +
  geom_smooth(method = "loess", se=T, color = "black",linetype = "dashed") +
  xlab("Years") + ylab("MAPP political party membership (in thousands)") +
  theme_light() +
  theme(legend.position="bottom") 

comparison.bw <- mappbw + opamed.plotbw
tiff("Figure 2 bw.tiff", units = "in", width = 10, height = 10, res = 350)
comparison.bw
dev.off()

cairo_ps("Figure 2 bw.eps")
comparison.bw
dev.off()

# Figure 5. Density plot ruling party size #### 
ridges1 <- opamed %>%
  mutate(party = fct_relevel(party, 
                             levels = "COLORADO", "PSUV", "CMM", "BTHSYR", "WPK", "PRI",
                             "KSC", "MPLA", "CPR", "BTH", "SED", "UMNO", "FRELIMO", "KMT", 
                             "AKP", "MDP", "BKP", "MKP", "LCY", "PSPR", "EPRFD", "CPCu",  
                             "GEORCP", "BLRCP", "CPSU", "CPLa", "CPEs", "CCP", "CPLi", "CPV", "AZERCP", "ARMCP",
                             "KZCP", "PLA", "MPRP","RPL", "MLDVCP", "UZCP", "KIRGCP", "TURKCP", 
                             "TADZCP", "UNITRUS", "PAP"))

fig.5 <- ggplot(ridges1) +
  geom_density_ridges(data = ridges1, aes(x = mempop_100, y = factor(party), 
                                          fill = factor(continent)), 
                      jittered_points = T, position = "raincloud",
                      alpha = 0.7, scale = 0.9) +
  xlab("Size") + ylab("Name") +
  scale_fill_manual(values = wes_palette("Darjeeling1", n = 5),
                    labels = c("Africa", "Latin America", "Asia", "Europe", "MENA")) +
  scale_x_continuous(breaks=seq(0, 40, by = 5)) +
  labs(fill = "Continent") +
  theme_ridges() + 
  theme(legend.position = "bottom")

tiff("Figure 5.tiff", units = "in", width = 10, height = 10, res = 350)
fig.5
dev.off()

cairo_ps("Figure 5.eps")
fig.5
dev.off()

# fig.5 black and white 
fig.5bw <- ggplot(ridges1) +
  geom_density_ridges(data = ridges1, aes(x = mempop_100, y = factor(party), 
                                          fill = factor(continent)), 
                      jittered_points = T, position = "raincloud",
                      alpha = 0.7, scale = 0.9) +
  xlab("Size") + ylab("Name") +
  scale_fill_grey(labels = c("Africa", "Latin America", "Asia", "Europe", "MENA")) +
  scale_x_continuous(breaks=seq(0, 40, by = 5)) +
  labs(fill = "Continent") +
  theme_ridges() + 
  theme(legend.position = "bottom")

tiff("Figure 5 bw.tiff", units = "in", width = 10, height = 10, res = 350)
fig.5bw
dev.off()

cairo_ps("Figure 5 bw.eps")
fig.5bw
dev.off()

# Figure 6 hegemonic and single-party size ####

regime.penetration <- ggplot(opamed, aes(regime.dummy, mempop_100, fill = factor(regime.dummy))) +
  geom_boxplot()+
  xlab("Regime") + ylab("Party size") +
  scale_fill_manual(values = wes_palette("Darjeeling1", n = 2),
                    labels = c("Hegemonic party", "Single party")) +
  scale_x_continuous(breaks=seq(0, 2, by = 1)) +
  labs(fill = "Regime division") +
  theme_light() +
  theme(legend.position="bottom") 

tiff("Figure 6.tiff", units = "in", width = 8, height = 8, res = 350)
regime.penetration
dev.off()

cairo_ps("Figure 6.eps")
regime.penetration
dev.off()

# Figure 6 bw 
regime.penetration.bw <- ggplot(opamed, aes(regime.dummy, mempop_100, fill = factor(regime.dummy))) +
  geom_boxplot()+
  xlab("Regime") + ylab("Party size") +
  scale_fill_grey(labels = c("Hegemonic party", "Single party")) +
  scale_x_continuous(breaks=seq(0, 2, by = 1)) +
  labs(fill = "Regime division") +
  theme_light() +
  theme(legend.position="bottom") 

tiff("Figure 6 bw.tiff", units = "in", width = 8, height = 8, res = 350)
regime.penetration.bw
dev.off()

cairo_ps("Figure 6 bw.eps")
regime.penetration.bw
dev.off()

# Figure 7. Membership volatility 1945-2020 ####
memb.growth.plot <- ggplot(subset(opamed, memb_growth < 50 & memb_growth > -20), 
       aes(x = year, y = memb_growth, color = factor(survived))) + 
  geom_line() +
  xlab("Years") + ylab("Membership stability") +
  geom_hline(yintercept = 0, color = "black", linetype = "dashed") +
  scale_color_manual(values = wes_palette("Darjeeling1", n = 3),
                     labels = c("Failed", "Survived regime change", "Survived")) +
  labs(color = "Ruling party & regime change") +
  theme_light() +
  theme(axis.text.x = element_text(angle = 90),
        legend.position="bottom") +
  facet_wrap(~party, scale="fixed")

tiff("Figure 7.tiff", units = "in", width = 12, height = 10, res = 350)
memb.growth.plot
dev.off()

cairo_ps("Figure 7.eps")
memb.growth.plot
dev.off()

# Fig 7 bw
memb.growth.plot.bw <- ggplot(subset(opamed, memb_growth < 50 & memb_growth > -20), 
           aes(x = year, y = memb_growth, linetype = factor(survived))) + 
  geom_line() +
  scale_linetype_manual(values = c(1:3),
                        labels = c("failed", "survived regime change", "survive"))+
  scale_color_manual(values = c("black", "black", "black"),
                     labels = c("failed", "survived regime change", "survive"))+
  xlab("Years") + ylab("Membership stability") +
  geom_hline(yintercept = 0, color = "black", linetype = "dashed") +
  labs(linetype = "Ruling party & regime change") +
  theme_light() +
  theme(axis.text.x = element_text(angle = 90),
        legend.position="bottom") +
  facet_wrap(~party, scale="fixed")

tiff("Figure 7 bw.tiff", units = "in", width = 12, height = 10, res = 350)
memb.growth.plot.bw
dev.off()

cairo_ps("Figure 7 bw.eps")
memb.growth.plot.bw
dev.off()

# Figure 8. relationship size and growth ####
rel.siz.grt <- ggplot(subset(opamed, memb_growth<150 & mempop_100<10), 
       aes(x =mempop_100, y =  memb_growth)) + 
  geom_point(size = .5) +
  geom_smooth(method='lm') +
  scale_y_continuous(breaks=seq(-80, 150, by = 5)) +
  scale_x_continuous(breaks=seq(0, 10, by = 1)) +
  stat_cor(label.x = 1.5, label.y = -75)+
  xlab("Ruling party size") + ylab("Membership volatility") +
  geom_hline(yintercept = 0, color = "black", linetype = "dashed") +
  theme_light() +
  theme(legend.position="bottom")

tiff("Figure 8.tiff", units = "in", width = 10, height = 8, res = 350)
rel.siz.grt
dev.off()

cairo_ps("Figure 8.eps")
rel.siz.grt
dev.off()

rel.siz.grt.bw <- ggplot(subset(opamed, memb_growth<150 & mempop_100<10), 
                      aes(x =mempop_100, y =  memb_growth)) + 
  geom_point(size = .5) +
  geom_smooth(method='lm',  color = "darkgrey", linetype = "dashed") +
  scale_y_continuous(breaks=seq(-80, 150, by = 5)) +
  scale_x_continuous(breaks=seq(0, 10, by = 1)) +
  stat_cor(label.x = 1.5, label.y = -75)+
  xlab("Ruling party size") + ylab("Membership volatility") +
  geom_hline(yintercept = 0, color = "black", linetype = "dashed") +
  theme_light() +
  theme(legend.position="bottom")

tiff("Figure 8 bw.tiff", units = "in", width = 10, height = 8, res = 350)
rel.siz.grt.bw
dev.off()

cairo_ps("Figure 8 bw.eps")
rel.siz.grt.bw
dev.off()

# Figure 9. History ruling parties' growwth ####
historical.1 <- ggplot(subset(opamed, foundation_miller<100), 
                       aes(x = foundation_miller, y = memb_growth)) + 
  stat_summary(geom = 'point') +
  stat_summary(geom = "line") +
  ylab("Membership volatility") +
  geom_hline(yintercept = 0, color = "black", linetype = "dashed") +
  theme_light()+
  theme(axis.title.x=element_blank(),
        axis.text.x=element_blank(),
        legend.position = "bottom")

historical.2 <- ggplot(subset(opamed, foundation_miller<100), 
                       aes(x = foundation_miller, y = mempop_100)) + 
  stat_summary(geom = 'point') +
  stat_summary(geom = "line") +
  xlab("Years since political party established") + ylab("Ruling party size") +
  geom_hline(yintercept = 0, color = "black", linetype = "dashed") +
  theme_light()+
  theme(legend.position = "bottom") 

historical.patch <- historical.1 / historical.2

tiff("Figure 9.tiff", units = "in", width = 10, height = 10, res = 350)
historical.patch
dev.off()

cairo_ps("Figure 9.eps")
historical.patch
dev.off()


# Figure 10. OPAMED strategy/dimension depending on ideology ####
marxist <- ggplot(subset(opamed, mempop_100 < 15 & memb_growth<200), 
       aes(x = mempop_100, y = memb_growth, color = factor(marxist))) + 
  geom_point(size = .6) +
  geom_smooth(method = "loess", se=T)+
  xlab("Ruling party size") + ylab("Membership volatility") +
  geom_hline(yintercept = 0, color = "black", linetype = "dashed") +
  scale_color_manual(values = wes_palette("BottleRocket2", n = 2),
                     labels = c("Non-Marxist", "Marxist")) +
  scale_x_continuous(breaks=seq(0, 15, by = 1)) +
  scale_y_continuous(breaks=seq(-80, 150, by = 10)) +
  labs(color = "Ideology") +
  theme_light()+
  theme(legend.position = "bottom") 

tiff("Figure 10.tiff", units = "in", width = 10, height = 8, res = 350)
marxist
dev.off()

cairo_ps("Figure 10.eps")
marxist
dev.off()

# fig. 10 bw
marxist.bw <- ggplot(subset(opamed, mempop_100 < 15 & memb_growth<200), 
                  aes(x = mempop_100, y = memb_growth, linetype = factor(marxist))) + 
  geom_point(size = .6, alpha = .5) +
  geom_smooth(method = "loess", se=T, color = "black")+
  xlab("Ruling party size") + ylab("Membership growth") +
  geom_hline(yintercept = 0, color = "black", linetype = "dashed") +
  scale_color_manual(values = c("black", "black"),
                     labels = c("Non-Marxist", "Marxist"))+
  scale_linetype_manual(values = c(1,2),
                     labels = c("Non-Marxist", "Marxist")) +
  scale_x_continuous(breaks=seq(0, 15, by = 1)) +
  scale_y_continuous(breaks=seq(-80, 150, by = 10)) +
  labs(linetype = "Ideology") +
  theme_light()+
  theme(legend.position = "bottom") 

tiff("Figure 10 bw.tiff", units = "in", width = 10, height = 8, res = 350)
marxist.bw
dev.off()

cairo_ps("Figure 10 bw.eps")
marxist.bw
dev.off()

# clear
rm(list=ls())
